/*******************************************************************************

***** Why Working from Home Varies Across Countries and People

* The following code produces the main figures

*******************************************************************************/


** Define path:
gl path "C:\Users\pz0511\Dropbox\GlobalWFH\WFH and Culture\PNAS Brief Report\Replication package"

use "${path}/GSWA_w4.dta", clear



** First, we collapse the data to get country-specific averages

collapse n_work_home1 idv l_gdp swaa_av_ind dens lsi tot_commuting_time , by(country order_countries)



***************************************************
** Figure 1: More individualistic countries work from home more
***************************************************

	* Define clock, for label position
		 cap drop clock
		 gen clock = 3
		 replace clock = 9 if inlist(country,"Argentina","Norway","Sweden","Ireland","Czech Rep.")
		 replace clock = 6 if inlist(country,"Mexico")
		 replace clock = 12 if inlist(country,"Netherlands","Singapore","Thailand")
		 
		 
		 twoway (scatter n_work_home1 idv, mlab(country) mlabvpos(clock)  msymbol(d)) (lfit n_work_home1 idv) , legend(off) ylabel(0(0.5)2) xtitle("Individualism") ytitle("Average Full Paid Days WFH") 

		 
		 graph export "${path}/PNAS_Figure1.png", as(png) replace
		 graph export "${path}/PNAS_Figure1.eps", as(eps) replace





***************************************************
** Figure 2: Cultural Individualism is the most important explanator of cross-country variation in WFH
***************************************************

** To run the regressions, we standardize the regressors to have mean of 0 and SD of 1:

gl XVARS idv l_gdp swaa_av_ind lsi dens tot_commuting_time

foreach var of global XVARS {
cap drop temp
egen temp = std(`var')
replace `var' = temp
}




sum $XVARS



* Main multivariate regression
reg n_work_home1 $XVARS , r
matlist r(table)


foreach j of numlist 2/7 {
	local k = `j' - 1
matrix Reg_values_`j' = J(3,7,.)
matrix Reg_values_`j'[1, `j'] = r(table)[1,`k']
matrix Reg_values_`j'[2, `j'] = r(table)[5,`k']
matrix Reg_values_`j'[3, `j'] = r(table)[6,`k']
}

matrix R2_1 = J(1,7,.)
matrix R2_2 = J(1,7,.)
matrix R2_3 = J(1,7,.)
matrix R2_4 = J(1,7,.)
matrix R2_5 = J(1,7,.)
matrix R2_6 = J(1,7,.)
matrix R2_7 = J(1,7,.)


matrix R2_1[1,1] = 100*e(r2)

matlist R2_1

matrix empty = J(3,7,.) // for visual purposes




*** Univariate regressions:
local j = 2
foreach var of global XVARS {
	reg n_work_home1 `var'
	matrix R2_`j'[1, `j'] = 100*e(r2)
	local j = `j' + 1
}








***************************************************
		 ** Panel A: Coefficient estimates
***************************************************
coefplot matrix(empty) matrix(Reg_values_2) matrix(Reg_values_3) matrix(Reg_values_4) matrix(Reg_values_5) matrix(Reg_values_6) matrix(Reg_values_7), ciopts(recast(rcap)) recast(bar) barwidth(0.7) bcolor(%70) ci((2 3)) mlabel(string(@b, "%5.3f")) mlabcolor(black) mlabsize(medsmall) offset(0) vertical coeflabels(c1 = " " c2 = "{bf:Individualism}" c3 = `" "{bf:GDP}" "{bf:Per Capita}"  "'  c4 = `" "{bf:Industry}" "{bf:Mix}"  "' c5 = `"  "{bf:Lockdown}"  {bf:Stringency}"' c6 = "{bf:Density}" c7 = `" "{bf:Commute}"  "{bf:Time}"  "') legend(off )  groups(c2 c5 = " ", gap(0.6))  order(c1 c2 c6 c4 c3 c7 c5) xsize(8) ysize(4) ytitle("Coefficient (Standardized)") yline(0)

		 		 graph export "${path}/PNAS_Figure2A.png", as(png) replace
		 		 graph export "${path}/PNAS_Figure2A.eps", as(eps) replace


***************************************************
		 ** Panel B: Variance explained
***************************************************
coefplot  matrix(R2_1) matrix(R2_2) matrix(R2_3) matrix(R2_4) matrix(R2_5) matrix(R2_6)  matrix(R2_7) , ciopts(recast(rcap)) recast(bar) barwidth(0.7) bcolor(%70) mlabel(string(@b, "%5.1f"))  mlabsize(medsmall) mlabcolor(black)   ///	
		  vertical coeflabels(c1 = "{bf:All variables}" c2 = "{bf:Individualism}" c3 = `" "{bf:GDP}" "{bf:Per Capita}"  "'  c4 = `" "{bf:Industry}" "{bf:Mix}"  "' c5 = `"  "{bf:Lockdown}"  {bf:Stringency}"' c6 = "{bf:Density}" c7 = `" "{bf:Commute}"  "{bf:Time}"  "' ) legend(off )  groups(c2 c5 = " ", gap(0.6))  order(c1 c2 c6 c4 c3 c7 c5  ) xsize(8) ysize(4) ytitle("Explained Variance (R2)") offset(0)

		 graph export "${path}/PNAS_Figure2B.png", as(png) replace
		 graph export "${path}/PNAS_Figure2B.eps", as(eps) replace
		 
		 
		 

		
		